Votre projet RStudio doit être structuré de la manière suivante :
Vous devez enregistrer votre script régulièrement dans le fichier exo1.R.
Source du fichier Occitanie.gpkg : ADMIN EXPRESS COG édition 2019, IGN
Source du fichier base_cc_comparateur.xls : Base comparateur de territoires 2019, INSEE
sf et la fonction st_read().st_crs(). S’agit-il de données projetée ? Si ce n’est pas le cas, transformez la couche des communes dans la projection française (Lambert 93, EPSG : 2154) avec la fonction st_transform().## Linking to GEOS 3.7.1, GDAL 2.4.0, PROJ 5.2.0
## Reading layer `Occitanie' from data source `/home/tim/Documents/prz/carto_avec_r_exo/data/Occitanie.gpkg' using driver `GPKG'
## Simple feature collection with 4454 features and 11 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -0.3271723 ymin: 42.33292 xmax: 4.845565 ymax: 45.04557
## epsg (SRID): 4326
## proj4string: +proj=longlat +datum=WGS84 +no_defs
## Coordinate Reference System:
## EPSG: 4326
## proj4string: "+proj=longlat +datum=WGS84 +no_defs"
st_geometry())Pour connaitre la liste de tous les noms ou code de région on peut utiliser la fonction unique().
## [1] "32" "31" "82" "12" "81" "30" "11" "65" "46" "09" "34" "66" "48"
st_union ().aggregate() pour regrouper les polygones et calculer les sommes des populations communales.reg76 <- st_union(occ)
dep76 <- aggregate(occ["POPULATION"], by = list(occ$INSEE_DEP), sum)
plot(st_geometry(occ), col = "lightblue1", border = "white", lwd = .5)
plot(st_geometry(dep76), col = NA, border = "lightblue2", lwd = 1, add = TRUE)
plot(reg76, col = NA, border = "lightblue3", lwd = 2, add = TRUE)st_buffer(). Quel est le code INSEE de la commune de Toulouse?toulouse <- com31[com31$INSEE_COM == "31555",]
toulouse_buffer <- st_buffer(toulouse, 10000)
plot(st_geometry(com31), lwd = .5)
plot(st_geometry(toulouse), col = "darkblue", add = TRUE)
plot(st_geometry(toulouse_buffer), col = NA, lwd = 3, border = "red",
add = TRUE)Déterminez quelles communes de la Haute-Garonne intersectent le buffer créé.
st_intersects().## Simple feature collection with 2 features and 12 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 504779 ymin: 6198834 xmax: 565151 ymax: 6270449
## epsg (SRID): 2154
## proj4string: +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
## ID STATUT INSEE_COM INSEE_ARR INSEE_DEP
## 4 COMMUNE_0000000009762857 Commune simple 31046 2 31
## 5 COMMUNE_0000000009760450 Commune simple 31547 1 31
## INSEE_REG CODE_EPCI NOM_COM_M POPULATION TYPE NOM_COM
## 4 76 200072635 BAREN 12 COM Baren
## 5 76 200068641 SEYSSES 8787 COM Seysses
## geom in_buffer
## 4 MULTIPOLYGON (((504930 6199... FALSE
## 5 MULTIPOLYGON (((556904 6267... TRUE
Affichez/superposez toutes les couches géographiques créees :
Jouez sur les styles pour les différencier…
plot(reg76, col = NA, border = "grey50", lwd = 2, bg = "lightyellow")
plot(st_geometry(com31), col = "#aec8f2", border = "white", lwd = .5,
add = TRUE)
plot(st_geometry(com31[com31$in_buffer == TRUE, ]),col = "red", lwd = .5,
border = "white", add = TRUE)
plot(st_geometry(toulouse), col = "darkblue", border = "black", lwd = 1,
add = TRUE)
plot(st_geometry(toulouse_buffer), col = NA, border = "black", lwd = 2,
lty = 2, add = TRUE)
plot(st_geometry(dep76), col = NA, border = "grey50", add = TRUE)Créer un objet sf (points) contenant la localisation de la préfécture de région à Toulouse.
st_points(), puis un objet sfc avec avec st_sfc() et finalement un objet sf avec st_sf().# Sur OSM : https://www.openstreetmap.org/search?whereami=1&query=43.59825%2C1.45047#map=19/43.59825/1.45047
# Position de la préfecture 43.59825, 1.45047
library(sf)
pref_pt <- st_point(c(1.45047, 43.59825))
pref_sfc <- st_sfc(pref_pt, crs = (4326))
pref <- st_sf(name = "Préfecture", geometry = pref_sfc)
pref## Simple feature collection with 1 feature and 1 field
## geometry type: POINT
## dimension: XY
## bbox: xmin: 1.45047 ymin: 43.59825 xmax: 1.45047 ymax: 43.59825
## epsg (SRID): 4326
## proj4string: +proj=longlat +datum=WGS84 +no_defs
## name geometry
## 1 Préfecture POINT (1.45047 43.59825)
Calculez une matrice de distance entre la préfecture et les centroïdes des communes du département Haute-Garonne.
Ajouter cette distance dans une nouvelle colonne dist_pref dans com31. Pour cela, utiliser les fonction st_centroid() et st_distance().
N’oubliez pas de vérifier les projections utilisées avec st_crs(), au besoin modifiez les avec st_transform().
## Coordinate Reference System:
## EPSG: 4326
## proj4string: "+proj=longlat +datum=WGS84 +no_defs"
## Coordinate Reference System:
## EPSG: 2154
## proj4string: "+proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
## [1] FALSE
## [1] TRUE
# Centroides des communes
com31_centro <- st_centroid(st_geometry(com31))
com31$dist_pref <- st_distance(com31_centro, pref)
head(com31, 2)## Simple feature collection with 2 features and 13 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 504779 ymin: 6198834 xmax: 565151 ymax: 6270449
## epsg (SRID): 2154
## proj4string: +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
## ID STATUT INSEE_COM INSEE_ARR INSEE_DEP
## 4 COMMUNE_0000000009762857 Commune simple 31046 2 31
## 5 COMMUNE_0000000009760450 Commune simple 31547 1 31
## INSEE_REG CODE_EPCI NOM_COM_M POPULATION TYPE NOM_COM
## 4 76 200072635 BAREN 12 COM Baren
## 5 76 200068641 SEYSSES 8787 COM Seysses
## geom in_buffer dist_pref
## 4 MULTIPOLYGON (((504930 6199... FALSE 104829.70 [m]
## 5 MULTIPOLYGON (((556904 6267... TRUE 17211.48 [m]
Cartographiez distance de chaque commune du département à la préfecture avec la fonction plot().
# Avec un peu de paramétrage
plot(com31["dist_pref"], main = "Distance à la préfecture (en mètres)",
pal = hcl.colors(13,"Turku", rev = TRUE), border = NA,
key.pos = 1,key.width = .15, key.length = .75, graticule = TRUE,
reset = FALSE)
plot(st_geometry(pref), pch = 20, col = "red", add = TRUE )osmdata extrayez l’ensemble des restaurants présents dans le département.## Data (c) OpenStreetMap contributors, ODbL 1.0. https://www.openstreetmap.org/copyright
# Définition d'une bounding box
q <- opq(bbox=st_bbox(st_transform(com46,4326)))
# Extraction des restaurants
res <- add_osm_feature(opq = q, key = 'amenity', value = "restaurant")
res.sf <- osmdata_sf(res)
res.sf.pts <- res.sf$osm_points[!is.na(res.sf$osm_points$amenity),]
resto <- st_transform(res.sf.pts, st_crs(com46))
resto <- st_intersection(resto, st_geometry(com46))## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries
# Affichage des restaurants
plot(st_geometry(com46), col="darkseagreen3", border="darkseagreen4")
plot(st_geometry(resto), add=TRUE, pch=20, col = "#330A5FFF", cex = 0.5)# Compter les restaurants par communes
inter <- st_intersects(x = com46, y = resto)
com46$nresto <- sapply(inter, length)
# communes sans restos
com46noresto <- com46[com46$nresto==0, ]
plot(st_geometry(com46), col="darkseagreen3", border="darkseagreen4")
plot(st_geometry(com46noresto), col = 'red', add = TRUE)
plot(st_geometry(resto), add=TRUE, pch=20, col = "blue", cex = 1)# index du restaurant le plus proche
index <- st_nearest_feature(x = st_centroid(com46noresto),
y = resto)## Warning in st_centroid.sf(com46noresto): st_centroid assumes attributes are
## constant over geometries of x
# distance au plus proche
com46noresto$dresto <- st_distance(x = st_centroid(com46noresto),
y = resto[index, ],
by_element = TRUE)## Warning in st_centroid.sf(com46noresto): st_centroid assumes attributes are
## constant over geometries of x
# Affichage de la carte
plot(com46noresto['dresto'], reset = F,
main = "Distance au restaurant le plus proche")
plot(st_geometry(com46), col=NA, add= TRUE)
plot(st_geometry(resto), add=TRUE, pch=20, col = "red", cex = 1)Dans le même projet RStudio créez un deuxième script nommé exo2.R.
Votre projet RStudio doit maintenant être structuré de la manière suivante :
sf et la fonction st_read() pour importer les données.st_crs().readxl et la fonction read_excel() pour ouvrir la table de données correctement. Importer la table dans un objet nommé occ_df.## Reading layer `Occitanie' from data source `/home/tim/Documents/prz/carto_avec_r_exo/data/Occitanie.gpkg' using driver `GPKG'
## Simple feature collection with 4454 features and 11 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -0.3271723 ymin: 42.33292 xmax: 4.845565 ymax: 45.04557
## epsg (SRID): 4326
## proj4string: +proj=longlat +datum=WGS84 +no_defs
## Coordinate Reference System:
## EPSG: 4326
## proj4string: "+proj=longlat +datum=WGS84 +no_defs"
occ <- st_transform(x = occ_raw, crs = 2154)
library(readxl)
occ_df <- read_excel(path = "data/base_cc_comparateur.xls", sheet = 1, skip = 5)
head(occ_df, 2)## # A tibble: 2 x 36
## CODGEO LIBGEO REG DEP P16_POP P11_POP SUPERF NAIS1116 DECE1116 P16_MEN
## <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 01001 L'Abe… 84 01 767 780 16.0 41 25 306
## 2 01002 L'Abe… 84 01 243 234 9.15 21 7 101
## # … with 26 more variables: NAISD18 <dbl>, DECESD18 <dbl>, P16_LOG <dbl>,
## # P16_RP <dbl>, P16_RSECOCC <dbl>, P16_LOGVAC <dbl>, P16_RP_PROP <dbl>,
## # NBMENFISC16 <dbl>, PIMP16 <dbl>, MED16 <dbl>, TP6016 <dbl>,
## # P16_EMPLT <dbl>, P16_EMPLT_SAL <dbl>, P11_EMPLT <dbl>, P16_POP1564 <dbl>,
## # P16_CHOM1564 <dbl>, P16_ACT1564 <dbl>, ETTOT15 <dbl>, ETAZ15 <dbl>,
## # ETBE15 <dbl>, ETFZ15 <dbl>, ETGU15 <dbl>, ETGZ15 <dbl>, ETOQ15 <dbl>,
## # ETTEF115 <dbl>, ETTEFP1015 <dbl>
Créez un nouvel objet sf à partir d’une séléction par attribut`:
plot(st_geometry(...)).Pour connaitre la liste de tous les noms ou code de région voir la fonction unique()
## [1] "32" "31" "82" "12" "81" "30" "11" "65" "46" "09" "34" "66" "48"
merge(). Quel est l’identifiant commun entre la table géo et la table attributaire?## Simple feature collection with 2 features and 11 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 504779 ymin: 6198834 xmax: 565151 ymax: 6270449
## epsg (SRID): 2154
## proj4string: +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
## ID STATUT INSEE_COM INSEE_ARR INSEE_DEP
## 4 COMMUNE_0000000009762857 Commune simple 31046 2 31
## 5 COMMUNE_0000000009760450 Commune simple 31547 1 31
## INSEE_REG CODE_EPCI NOM_COM_M POPULATION TYPE NOM_COM
## 4 76 200072635 BAREN 12 COM Baren
## 5 76 200068641 SEYSSES 8787 COM Seysses
## geom
## 4 MULTIPOLYGON (((504930 6199...
## 5 MULTIPOLYGON (((556904 6267...
## # A tibble: 2 x 36
## CODGEO LIBGEO REG DEP P16_POP P11_POP SUPERF NAIS1116 DECE1116 P16_MEN
## <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 01001 L'Abe… 84 01 767 780 16.0 41 25 306
## 2 01002 L'Abe… 84 01 243 234 9.15 21 7 101
## # … with 26 more variables: NAISD18 <dbl>, DECESD18 <dbl>, P16_LOG <dbl>,
## # P16_RP <dbl>, P16_RSECOCC <dbl>, P16_LOGVAC <dbl>, P16_RP_PROP <dbl>,
## # NBMENFISC16 <dbl>, PIMP16 <dbl>, MED16 <dbl>, TP6016 <dbl>,
## # P16_EMPLT <dbl>, P16_EMPLT_SAL <dbl>, P11_EMPLT <dbl>, P16_POP1564 <dbl>,
## # P16_CHOM1564 <dbl>, P16_ACT1564 <dbl>, ETTOT15 <dbl>, ETAZ15 <dbl>,
## # ETBE15 <dbl>, ETFZ15 <dbl>, ETGU15 <dbl>, ETGZ15 <dbl>, ETOQ15 <dbl>,
## # ETTEF115 <dbl>, ETTEFP1015 <dbl>
## Simple feature collection with 2 features and 46 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: 526058 ymin: 6245429 xmax: 588475 ymax: 6256838
## epsg (SRID): 2154
## proj4string: +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
## INSEE_COM ID STATUT INSEE_ARR INSEE_DEP
## 1 31001 COMMUNE_0000000009761264 Commune simple 2 31
## 2 31002 COMMUNE_0000000009761234 Commune simple 3 31
## INSEE_REG CODE_EPCI NOM_COM_M POPULATION TYPE NOM_COM LIBGEO REG DEP P16_POP
## 1 76 200072643 AGASSAC 115 COM Agassac Agassac 76 31 115
## 2 76 200071298 AIGNES 250 COM Aignes Aignes 76 31 250
## P11_POP SUPERF NAIS1116 DECE1116 P16_MEN NAISD18 DECESD18 P16_LOG P16_RP
## 1 115 9.58 5 10 51.0000 0 1 71.0000 51.0000
## 2 238 21.81 11 9 101.4073 5 0 122.4073 101.4073
## P16_RSECOCC P16_LOGVAC P16_RP_PROP NBMENFISC16 PIMP16 MED16 TP6016
## 1 20 0 38.00000 42 NA 20809.33 NA
## 2 10 11 88.35487 105 NA 22740.67 NA
## P16_EMPLT P16_EMPLT_SAL P11_EMPLT P16_POP1564 P16_CHOM1564 P16_ACT1564
## 1 22.96553 8.965533 31.01992 80.000 8.000000 55.0000
## 2 27.07635 8.030424 28.99414 154.619 9.036294 120.4819
## ETTOT15 ETAZ15 ETBE15 ETFZ15 ETGU15 ETGZ15 ETOQ15 ETTEF115 ETTEFP1015
## 1 18 7 0 3 6 2 2 4 0
## 2 21 9 3 3 4 0 2 4 0
## geometry
## 1 MULTIPOLYGON (((526952 6255...
## 2 MULTIPOLYGON (((584462 6253...
Réalisez huit cartes thématiques finalisées.
Chacune des cartes devra comprendre les éléments d’habillage et de mise en page nécessaires leur compréhension (titre, légende, sources…)
Parmi ces huits cartes deux échelles différentes (échelle de la région, échelle d’un département) et deux découpages territoriaux différents (les communes, les départements) doivent être utilisés.
Vous devez réaliser les huit types de cartes suivants :
Exemple :
Nous vous conseillons d’utiliser le package cartography pour l’ensemble de vos réalisations. Néanmoins l’utilisation d’autres packages comme ggplot2 ou tmap est possible.
Vous devrez nous rendre pour le vendredi 28 février 2020 un script permettant de construire ces huit cartes dans un fichier dont le nom est constitué de votre prénom et de votre nom (Prenom_Nom.R).
Votre script doit pouvoir fonctionner une fois inséré dans un projet RStudio structuré de la manière suivante :
Vous serez évalués de la manière suivante :
Pour répondre à cet exercice vous devrez créer au moins deux nouvelles variables. Une variable quantitative relative :
Une variable qualitative (n’utilisez par cette variables dans votre exercice ) :
# calcul de la part des résidences secondaires
# dans les logements
com31$logsec <- 100 * com31$P16_RSECOCC / com31$P16_LOG
summary(com31$logsec)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 1.707 5.626 13.164 16.059 88.889
com31$logsectypo <- cut(com31$logsec, breaks = c(0, 1.8, 16, 89),
labels = c("faible", "moyenne", "forte"),
include.lowest = TRUE)
typoLayer(com31, var = "logsectypo", col= c("red", "yellow", "green"),
legend.pos = "topleft",
legend.title.txt = "Part des résidences secondaire")# créer une couche des départements
dep76 <- aggregate(occ["POPULATION"], by = list(INSEE_DEP = occ$INSEE_DEP), sum)
library(cartography)
library(mapinsetr)
# créez un masque à la dimension de la région
mask <- create_mask(bb = st_bbox(st_buffer(dep76, 10000)))
# déplacer et redimensioner la région
minidep76 <- move_and_resize(x = dep76, mask = mask,
xy = c(564922.8, 6184225), k = 0.1 )
# déplacer et redimensioner le masque
minimask <- move_and_resize(x = mask, mask = mask,
xy = c(564922.8, 6184225), k = 0.1 )
par(mar = c(0, 0, 1.2, 0))
plot(st_geometry(com31), bg = "cornsilk", col = "lightblue4",
border = "lightblue2")
plot(st_geometry(minimask), col = "cornsilk2", border = NA, add = T)
plot(st_geometry(minidep76), col = "grey70", border = NA, add=T)
plot(st_geometry(minidep76[minidep76$INSEE_DEP =="31", ]),
col = "lightblue4", border = NA, add=T)
layoutLayer("Le département Haute-Garonne",
author = "ADMIN EXPRESS COG édition 2019, IGN",
sources = "T. Giraud, 2020", scale = 10,
horiz = F, posscale = "bottomleft", south = TRUE)